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Based on computer simulations and theoretical analysis, a new dynamic state is found in induc- 
tively coupled intrinsic Josephson junctions in the absence of an external magnetic field. In this 
state, the plasma oscillation is uniform along the c axis with the fundamental frequency satisfying 
the ac Josephson relation. There are (2m + 1)tt phase kinks around the junction center, with m 
being an integer, periodic and thus non-uniform in the c direction. In the IV characteristics, the 
state manifests itself as current steps occurring at all cavity modes. Inside the current steps, the 
plasma oscillation becomes strong, which generates several harmonics in frequency spectra at a given 
voltage. The recent experiments on terahertz radiations from the mesa of a BSCCO single crystal 
can be explained in terms of this state. 

PACS numbers: 74.50.+r, 74.25.Gz, 85.25.Cp 



Soon after the discovery of high-T c superconductivity 
in cuprates, it was demonstrated that a single crystal 
of layered superconductors, such as BSCCO, behaves as 
a stack of intrinsic Josephson junctions (IJJs) [1]. The 
IJJs have been investigated extensively since then, and 
various dynamical properties were explored. One focus 
in this field is to exploit the dynamic states of IJJs to 
generate electromagnetic waves in the terahertz (THz) 
regime 0, H, 0] ■ The high-T c cuprate superconductors are 
unique for this purpose since the IJJs are homogeneous at 
atomic scale and the energy gap is large, 60meV, which 
make it possible to tune the frequency in a wide range by 
adjusting the bias voltage. To obtain coherent radiations, 
however, all the junctions should be operated in-phase, 
which turns out to be not easy. 

The dynamics of IJJs are appropriately described by 
the coupled sine-Gorden equations [5| 

9 2 P ; = (1 - CA^XsinP; + f3d t Pl + dfP, - Jext), (1) 

where P/ is the gauge-invariant phase difference at the 
Zth junction, (3 = 4ira c \ c / ' c^[e~ c the normalized c-axis 
conductivity, £ = \ 2 ab /sD the inductive coupling and 
A( 2 )p = p +1 + P;_ x - 2P; the difference operator. 
e c is the dielectric constant and o~ c is the conductance 
along the c axis. A c and X a b are the penetration depths. 
s (D) is the thickness of the superconducting (insulat- 
ing) layer. Length is normalized by A c , time is nor- 
malized by the Josephson plasma frequency l/2nfj with 
fj = c/27rA cv / e^, and voltage is in unit of $o/j/c Here 
c is the light velocity in vacuum. 

Without an external in-plane magnetic field, Eq. (Q]) 
have the plasma solution, moving soliton and anti-soliton 
states, and McCumber state [6(. The plasma solution 
corresponds to the collective small phase oscillation in 
Eq. It resonates in the cavity formed by IJJs when 
the cavity resonance condition meets, i.e. q\ w /2 = L 



with q being an integer, X w plasma wavelength and L the 
length of the IJJs. It corresponds to the voltage per layer 
V = qncj/L with Cj the plasma velocity (1 < j < N, N 
the total number of junctions) Q assuming that the ac 
Josephson relation is satisfied. 

In the soliton state, a soliton (anti-soliton) moving 
with velocity v is reflected as an anti-soliton (soliton) by 
the edges. The periodic motion and reflection of soliton 
and anti-soliton contributes to a dc voltage V = 2irv/L. 
When the velocity of soliton approaches Cj ,it_gives birth 
to the zero-field step (ZFS) at V = 2-KCj/L ^AWk, which 
corresponds to the 2nd cavity mode. If n solitons engage 
in this process, they produce the nth ZFS. The inter- 
play between solitons and plasma leads to many peculiar 
behaviors [loj ]. such as Cherenkov radiation [ill ]. 

In the McCumber state, the phase is uniform through 
the x direction. Therefore the junctions are decoupled 
and it gives linear branches in the IV characteristics. 

In the present study, we find that in the absence of 
an external magnetic field, a new dynamic state can be 
achieved in the inductively coupled sine-Gordon systems, 
which is characterized as follows: in the IV characteristics 
it generates current steps at all cavity modes; while the 
plasma oscillation is uniform along the c axis, in the total 
phase there are static kinks of (2m + 1)tt along the x axis, 
which are periodic but non-uniform in the c direction. It 
is shown that the recent experiments 0, 0| on terahertz 
emission can be understood in terms of this new state. 

The results shown below are for (3 — 0.02, ( = 4/9 x 
10 6 , typically for BSCCO, and L = 80/im; confirmations 
of the main results have been made with other param- 
eters. As displayed in the inset of Fig. [TJ the current 
J cxt is fed in along the c axis. Since we are searching 
for states periodic in the c axis, we adopt the periodic 
boundary condition along this direction which makes the 
maximum plasma velocity cq = 1 possible. The equa- 
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FIG. 1: (color). Simulated IV characteristics under the 
boundary conditions (a) without and (b) with radiation with 
C = 0.000177 and R = 707.1. The vertical dashed lines cor- 
respond to the voltage at cavity modes V = n$o/(2y / eTL). 
Inset: geometry of the mesa of BSCCO single crystal. 



tions of motion Eq. {1} are integrated by the staggered 
leapfrog algorithm. The time step in all simulations is set 
to At = 0.0018 and the mesh size is set to Ax = 0.002. 
The accuracy is checked with smaller Ax and At. We 
consider first the system under the boundary condition 
d x Pi — 0, which turns out to be a good unperturbed 
state for discussion on radiation as revealed later. 

The important feature of the IV characteristics in Fig. 
[T]Ja) is the existence of current steps. In sharp contrast 
with the ZFS in single junctions [8j, the present current 
steps in IJJs occur at both even and odd cavity modes; 
the dc voltage and the frequency of plasma oscillation 
satisfy the ac Josephson relation. 

As revealed by the present numerical simulations, there 
are four parts in the phase variables: (1) a term evolving 
linearly with time according to the ac Josephson relation 
Lot, (2) a term of (2m + l)7r phase kinks and anti-kinks ar- 
ranged periodically in the c direction, (3) a plasma term 
uniform along the c axis P(x,t), (4) a constant term ip. 
Figure [2] exhibits several typical configurations for the 
2nd term. For ease of discussions, we restrict ourselves 
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FIG. 2: (color). Typical configurations of the static term in 
the phase difference (the 2nd term in Eq. lHJ): (a) period of 
2 layers and (b) period of 4 layers, corresponding to the 1st 
current step, and (c) corresponding to the 2nd current step. 



to the 1st cavity mode in the following discussions, which 
can be extended to higher modes straightforwardly. The 
total phase therefore can be expressed by 



Pi(x,t) =ut + fiP l s (x)+P(x,t)+ip, 



(2) 



where P*(x) is a soliton from to 27r with the center at 
x = L/2 (see below), P(x, t) = A cos(irx / L) s'muit+ high 
harmonics. 

Substituting Eq. © into Eq. (Q]) and using the ex- 
pansion of sine of sine by the Besscl function, we arrive 
at the following relation for the static parts of phases 



fidlPf =U0- Jcxt 

-(1 - CA (2) ) Ji(Acos(7o;/L)) sin(/,P ; s + ip), 



(3) 



where J\ is the Bessel function of the first kind and for 
simplicity only the fundamental mode is taken. 

Off resonance the system behaves in an ohmic way as 
seen in Fig. [T] where J cxt ~ /3ui and ip w 0. The solutions 
of period of 2 layers (period-2), where /; = (—1)'/ (Fig. 
E£a)), and of 4 layers (period-4), where /; = (— 1)['/ 2 !/ 
(Fig. [2fb) ) , diagonalize the difference operator A^ 2 ). In 
both cases, Pf is independent of I, thus Pf = P s0 . We 
then arrive at the equation 

f0 2 x P s0 = -(l + q OJi(Acos(nx/L)) S in(fP s0 ), (4) 

where q — 4 for period-2 and q = 2 for period-4. 

Since the l.h.s. of Eq. (|4]) and J\(Acos(x-k / L)) are 
antisymmetric with respect to x — L/2, sin(/P s0 ) must 
be symmetric. This imposes / = (2m + l)/2, with an 
integer m. 

When the voltage is tuned to one of the cavity modes, 
the plasma oscillation can acquire large amplitude in the 
form of standing waves. By taking a non-zero constant 
tp, the system can bypass increasing J ox t through the 
Josephson current 

J cxt =(3uj- s -^ r Q L dx.h{Acos{irx/L)) cos(fP sQ ). 

(5) 

This gives rise to the current steps in the IV characteris- 
tics in Fig. [TJ Extra amount of energy besides the ohmic 
dissipation is fed into the plasma oscillation. 

In Figs. [3] (a)-(c), we display distributions of 
the dc part of Josephson current (sinP/) t in two 
adjacent junctions. The symbols are given by 
J 1 (Acos(irx/L))sin(±P s0 {x)/2 + ip) where P s0 (x) = 
4arctan[exp((a; - L/2)/\'j)] (see Eq. (@J|). The depen- 
dence of net Josephson current J nc t = (sinP;)^,* on J ex t 
at the 1st current step is shown in Fig.[3jd). 

The remaining terms in the r.h.s. of Eq. ([3]), which 
are proportional to the inductive coupling £, should cor- 
respond to a (unquantized) soliton as required by the l.h.s 
of the equation. This can be seen clearly in Fig. [3je). 
The system arranges itself in a way that the parameter 
C is blocked from involving directly in the net Josephson 
current, which should be in order of J c . 
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FIG. 3: (color). Comparison between the theory and simu- 
lations on the dc part of Josephson current (sinP;)t in two 
neighboring layers for the configuration in Fig. [2jb): (a) off 
resonance at J ex t = 0.12, (b) at the middle of the 1st step with 
Jcxt = 0.35, (c) at the top of the 1st step with J cxt = 0.53. 
The solid curves are results by numerical simulations and the 
symbols are derived by the analytic expressions, (d) Cur- 
rent dependence of voltage (red) and J nc t (blue) at the 1st 
current step. Green symbols are for the theoretical estimate 
J cxt — uif3. (e) {A^ sin Pi) t at J oxt = 0.35 (dashed red line) 
and at J cxt = 0.53 (dashed green line). 



Despite that the total phase is non-uniform along the 
c axis, the electric field Ef — d t Pi as well as the mag- 
netic field Bf = (1 - (A^)- 1 d x Pi = d x P + 0(1/0 
with Pi(x,t) given by Eq. ^ are uniform along the 
c axis. This makes a coherent superradiation possible as 
observed in the recent experiments [3J, |4j . Distributions 
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FIG. 4: (color). Time evolutions of the magnetic field and 
electric field and the frequency spectrum. The dc part of 
the electric field has been subtracted. Left/right: at the bot- 
tom/top of the 1st current step. 



of electric and magnetic fields are displayed in Fig. [4] In 
contrast with the state at bottom of the current step (left 
of Fig. 2]), the amplitude of plasma term increases and 
high harmonics become visible inside the current steps, 
and thus the standing waves of electromagnetic field are 
not simply trigonometric anymore (right of Fig. 2]). The 
wave crest of magnetic field becomes smaller when it ap- 
proaches the edges and it changes sign when it is re- 
flected. It resembles the motion of solitons in ZFS ex- 
cept that the magnetic flux is not quantized since it is 
associated with the plasma oscillations. 

From the foregoing minimum theoretical analysis, we 
have already captured the main features of the system 
revealed by simulations. The following discussions are 
necessary to make the picture more comprehensive: (I) 
There are infinite number of configurations of the static 
phase even in the same current step. Besides those in 
Figs. [2ta) and (b), periodic configurations of (2m + 1)tt 
kinks, and their composites are observed in simulations. 
Configurations of phase kinks not periodic in the c axis 
may be possible. Superpositions of soliton(s) and anti- 
soliton(s) along the junction direction is also permitted 
as for the single sine-Gordon system [13]. The height 
of current steps depends on the specific configuration as 
can be read from Eq. |(SJ). (II) The phase configura- 
tions for higher cavity modes can be constructed from 
the first cavity mode, see Fig. EJd). (Ill) In addition 
to the plasma oscillation, there exists a small oscillation 
of the center of soliton Pf(x) as well, which produces 
a spike-like oscillation with characteristic length X'j in 
electromagnetic fields. (IV) As revealed in the above dis- 
cussions, the strong inductive coupling is essential for 
the new dynamic states, and therefore it is impossible 
to build such current steps in a single Josephson junc- 
tion. In fact, the current steps completely vanish when 
£ is small in simulations. (V) While BSCCO may be the 
best choice, the new state can be realized in other layered 
superconductors, including s-wave ones with not so high 
T c . The crystal is to be thick to reduce the surface effect. 

Let us investigate the possible radiation of energy from 
the IJJs. The boundary condition should be formulated 
explicitly for this purpose. Theoretical efforts have been 
attempted to model the interaction between IJJs and 
the outside space [l3|, 14, 15, [l(|. We model the ex- 

and hence 

e*/b» = 



14] 



ternal load by an effective RC circuit 
the boundary condition becomes dtP/d x P 
Z = R — i/(Cuj) = \Z\exp(i8). The sharp radiation 
peaks upon voltage sweeping observed in recent experi- 
ments 0, IH indicate that there is a significant mismatch 
in impedance between IJJs and the outside space, and 
that electromagnetic waves are reflected mostly at the 
edges of IJJs. The small ratio between the mesa thick- 
ness in the experiments and A c may be a source of this 
mismatch [l6j], which implies \Z\ 3> 1. We have con- 
firmed that the new states remain stable for large values 
of \Z\. 
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FIG. 5: (color). Dependence of radiation power evaluated by 
the Poynting vector at the edges of IJJs (a) on \Z\ and (b) 
on 6 inside the 1st current step at J oxt = 0.27. (c) Radiation 
power at all current steps. 

The IV characteristics under the boundary condition 
permitting radiation are shown in Fig. QJb). The heights 
of the current steps are reduced when radiations are 
present, since when the radiation is strong it is hard to 
form stable standing waves. 

As shown in Fig. EIc), the radiation power increases 
linearly with the external current at the current steps. 
This behavior can be explained by the above theoreti- 
cal description including additionally the radiation as a 
small perturbation. The power at the top of the 1st cur- 
rent step is of one order of magnitude stronger than that 
at the bottom, and the maximum power with the present 
boundary condition is about 8000 W/cm 2 , which is quite 
promising for applications. The dependence of the radia- 
tion power on \Z\ and 9 is detailed in Figs. [UJa) and (b), 
which can be described as cos9/\Z\. 

The above simulation results seem to be able to pro- 
vide a consistent picture for the two recent experiments 
0, Q] on the terahertz radiation in zero external magnetic 
field. Coherent radiations were detected in the resistive 
curve in Ref. Q with the frequency corresponding to the 
fundamental mode of the 1st current step. In Ref. [3], 
radiations were detected in the anomalous region of IV 
curve, and the detected power is much stronger than that 
in Ref. [3( . Comparing these features with our simulated 
results, it is natural to assign the radiations in Ref. [3| to 
states inside the 1st current step, while those in Ref. Q 
to states at the bottom of the current step. The shape 
of the anomalous IV curve where radiations are observed 
in Ref. 0] is also indicative that the system falls into the 
current step. High cavity modes could not be observed 
in experiments due to the heating effect, which is not 
taken into account in the present model. The temper- 
ature dependence of the conductivity j3 is necessary in 



order to simulate the overall shape of IV characteristics 
in the experiments especially at high currents. 

In contrast, the plasma solution without static kinks 
(crossover to the McCumber state without radiation) 
cannot show the sharp peak in radiation energy when 
voltage is tuned, and there is no current step because the 
plasma do not couple to the cavity modes. On the other 
hand, the radiation from the soliton state does not sat- 
isfy the ac Josephson relation, and thus it cannot explain 
the recent experiments. 

In conclusion, based on computer simulations and the- 
oretical analysis we reveal the existence of a new dy- 
namic state in inductively coupled intrinsic Josephson 
junctions. In this state, the plasma oscillation is uniform 
along the c axis with the fundamental frequency satisfy- 
ing the ac Josephson relation, while static phase kinks of 
(2m + l)7r are stacked periodically. The states manifest 
themselves in the IV characteristics as current steps at 
all cavity modes. This state supports quite strong radi- 
ations. The recent experiments on THz radiations from 
BSCCO single crystals can be interpreted in terms of this 
state in a consistent way . 
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